% by clq ylx
% 2021/6/12
%function [best_is_pref,profits_PP]=main_pref(P,LAMBDA,ci,cI,r,J,I,j_all)
function [best_is_pref,profits_PP]=main_pref(j_all)
    global ci cI r J I LAMBDA P
    addpath functions;
    Line = j_all;
    patient_number = sum(Line > 0); % equivalent to size(2)
    genPoisson(patient_number);
    genPIE(patient_number);
    for ii=1:100
        disp(ii)
        S_pref = zeros(I, J);
        % Profit = zeros(patient_number, I); % For What?

        init_time = ftimer(); % set main clock

        best_is_pref = zeros(1, patient_number);
        profits_pref = zeros(1,patient_number+1);
        pref=randi([1,2],1,patient_number);
        for n = 1:patient_number
            
            j = Line(n);
            profits = zeros(1, I);

            for i = 1:I
                S_test = S_pref;
                S_test(i, j) = S_test(i, j) + 1; % Replace rim & rimS
                Q = findQ(S_test); 
                R = findR(Q);
                profits(i) = function_fQR(Q, R);
            end
            
        %% special part for pref
            pref_temp=pref(n);
            if pref_temp==1
                profit_temp=max(profits(1:I/2));
                i=find(profits == profit_temp);
            else
                profit_temp=max(profits(I/2+1:end));
                i=find(profits == profit_temp);
            end
            profits_pref(n+1) = profit_temp + profits_pref(n);
            best_is_pref(n) = i;
            S_pref(i, j) = S_pref(i, j) + 1;
            Q = findQ(S_pref);
            R = findR(Q);
            profits_PP(ii,n) = profits(i);  %% PP refers to proposed preference
            % fprintf("    Time = %.2f seconds \n\n", ftimer()-start_time);
        end
    end
    
    profits_PP = mean(profits_PP,1);
            

        fprintf("\n\n Totle Time = %.1f seconds \n", ftimer()-init_time);

        display(S_pref)
        display(best_is_pref)

        % plot(1:patient_number,sum_profits);
        %showTimer();